Bunching Transitions on Vicinal Surfaces and Quantum N-mers 
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We study vicinal crystal surfaces with the terrace-step-kink model on a discrete lattice. Including 
both a short-ranged attractive interaction and a long-ranged repulsive interaction arising from elastic 
forces, we discover a series of phases in which steps coalesce into bunches of n steps each. The value 
of n varies with temperature and the ratio of short to long range interaction strengths. We propose 
that the bunch phases have been observed in very recent experiments on Si surfaces. Within the 
context of a mapping of the model to a system of bosons on a ID lattice, the bunch phases appear 
as quantum n-mers. 
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The study of the equilibrium properties of a stepped 
vicinal crystal surface is important from technological 
and fundamental perspectives. It also provides a fasci- 
nating example of a problem which has a much broader 
context, reaching, via a mapping onto a one-dimensional 
quantum chain, into the realm of ID quantum liquids. In 
this Letter we address two related crystal surface prob- 
lems. The first is that of the description of apparent 
tricritical phenomena observed in a beautiful set of ex- 
periments by Song et al (!]] on silicon surfaces, miscut 
away from the [113] direction towards [001]. The second 
is the related very recent observations of multiple-height 
steps by Sudoh et al J2| on vicinal silicon surfaces near 
the [113] crystalline direction. 

Recent theoretical work [||J] attempted to explain the 
results of Song et al in terms of a continuum model of 
steps interacting via a long-ranged repulsive elastic in- 
teraction and a short-ranged attractive interaction. This 
model fails to describe the observed bunching of steps 
on the vicinal surfaces coexisting with the (113) facet. 
Here we explore the consequences of retaining the dis- 
crete, atomic nature of steps on a crystal surface within 
the same model. We discover entirely different physics. 
The steps do not phase separate but instead can coalesce 
into bunches whose size depends on the relative strengths 
of the short- and long-range interactions. Vicinal sur- 
face phases can be characterized by widely separated n- 
bunches and transitions occur between phases having dif- 
ferent values of n. We propose that bunch phases with 
large n correspond to the two-phase coexistence region 
of Song et al and that the n=2,3,4 bunch phases produce 
the multiple-height steps seen by Sudoh et al. 

We employ the terrace step kink model, which is a 
lattice model that can be mapped onto an equivalent 
ID quantum mechanical model of interacting spinless 
fermions or hard-core bosons ||. In the quantum pic- 
ture it is well-known that the 1-bunches form a Luttinger 



lattice liquid. Our results generalize this picture to a 
Luttinger lattice liquid which can form dimers (the 2- 
bunches), trimers (the 3-bunches), and in general n-mers. 
The transitions between these phases are quantum many- 
body phenomena of a type hitherto unexplored. 



Consider a crystal with lattice constants (a x 



where the steps are in the x — y plane and run, on the 
average, in the y— direction. In terms of hard-core boson 
operators our Hamiltonian, for a system with L lattice 
sites in the x— direction, is written as 
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where t, a monotonically increasing function of the tem- 
perature T, is a hopping matrix element related to the 
step stiffness £ by £ = - „ B 2 1 , rj{T) is the free energy 



per unit length of the steps, G measures the strength of 
the elastic interaction Q between steps, and U gives a 
short range attraction between steps. Here a^'s and a\ 
are hard-core bosonic creation and annihilation operators 
respectively and n, = a\ai is the boson number operator 
on site i. This representation is equivalent to the fermion 
representation introduced in Ref. J^], as demonstrated 
by application of the unitary Jordan- Wigner transfor- 
mation M. For N bosons or steps confined on L-sites 
the ground state energy of H, E (s,T), is related to 
the surface energy of the vicinal surface, 7(s), through 
r j{s)/cos{9) = 70 + E (s,T)/La x , where s = N/L is the 
step or boson density, = tan _1 (s) is the miscut angle 
and 7q is the surface tension of the reference surface. 

We now show that a number of features of this many 
body system can be inferred from exact diagonalization 
of small systems in conjunction with some simple analyt- 
ical calculations. We find that for sufficiently attractive 



1 



interactions, the steps on the surface rearrange them- 
selves into bunches at low temperatures. At T = 0, 
where entropic effects are unimportant, the bunch size 
is completely determined by the ratio of the strengths of 
the attractive and repulsive interactions i.e. U/G. With 
increasing temperature, the steps start peeling off from 
the bunches in a series of bunching transitions, until only 
one step is left in a bunch at high temperatures. In what 
we call the extreme dilute limit, we obtain a phase dia- 
gram that shows these bunching transitions as function 
of temperature. 

Consider first the limit t — > 0, where the hopping part 
of the Hamiltonian Eq. (|l|) can be completely ignored; 
this is the zero temperature limit of the vicinal surface. 
In this limit, the energy of the system is obtained by 
minimizing the total potential energy, and the minimum 
energy configuration consists of bosons in bunches of size 
n&, well separated from each other. In a bunch, the 
bosons sit next to one another. For a bunch of size rib, 
the energy per boson or free energy per step is given by 

tl n s , n s U(n b -1) G \^ (n b - i) 

f(n b , 0) = 17(0) + — > - 2 , 2 

m rib i 

where we have neglected contributions from bunch-bunch 
interactions that are smaller than the leading order terms 
given in Eq. (Q) by a factor s 2 <C 1. We call the limit 
in which these interactions are completely ignored, the 
"extreme dilute limit "(EDL). With increasing U/G, the 
bunch size keeps increasing; for U/G < 1 the bunch 
size is 1, for 1 < U/G < 1.5 the bunch size is 2, for 
1.5 < U/G < 1.83 the bunch size is 3, for 1.83 < 
U/G < 2.08 the bunch size is 4, and so on. At the 
points U/G = 1,1.5,1.83,2.08,... there is a coexistence 
of bunches that differ in size by one. From plots of the 
energy per boson vs. U/G for different rib's, we find that 
the slopes of these curves at the points where they inter- 
sect are not the same. The T = bunching transitions 
brought about by changing U/G can then be considered 
first order. 

At finite temperatures we expect the free energy per 
step in larger bunch sizes to decrease less rapidly com- 
pared to smaller bunches, as steps within larger bunches 
will have their fluctuations more constrained. As a result, 
larger bunches should rearrange themselves into smaller 
bunches, with increasing temperature, until eventually 
at very large temperatures only one step is left in the 
bunch. In order to study these bunching transitions, we 
introduce g — G/t and u = U/t, which are the interaction 
strengths scaled by the hopping matrix element t. Note 
that <7 _1 plays the role of temperature. As in the T = 
case, we look at the free energy per step in the n^-bunch 
phase, f(rib,T), obtained by exactly diagonalizing H for 
t 0. In order to obtain the phase boundaries between 
the rib-bunch phase and the n b — 1 phase for a fixed u/g, 
one has to find the value of g that satisfies the equation 



f(nb,T) — f(rib — l,T) = 0. We do this numerically by us- 
ing a standard root finding algorithm which approaches 
the root iteratively. The root finding algorithm requires 
repeated evaluation of f(n b ,T) — f(n b — 1,T) using ex- 
act diagonalization of the rib- and (n b — l)-bunches for 
different values of g. The results are shown in Fig. [I]. 
A key point emerging from our calculations is that, at 
the phase boundary between the n b - and (rib — l)-bunch 
phases, d[f(n b , T) - f(n b - 1, T)]/dT ^ 0, implying that 
the bunching transitions are again first order in the EDL. 

Fig, 1 
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FIG. 1. Phase diagram in the extreme dilute limit plotted 
in u/g-l/g space. Here 1/g plays the role of temperature. For 
a given surface, the ratio of the strengths of the short-range 
and long-range interactions fixes u/g, and increasing the tem- 
perature corresponds to increasing 1/g. The lines separate 
regions with stable bunch sizes that differ by one as indicated 
in the figure. The points marked at u/g = 2 correspond 
to the values of g for which the pair correlation functions, 
computed using Green's function Monte Carlo (GFMC), are 
shown in Fig. ^| For each value of g we also indicate the bunch 
size (lb, 2b, 3b and 4b correspond to l-,2-,3- and 4-bunches re- 
spectively, while cr refers to the critical region) obtained from 
GFMC simulations. 

In performing the exact diagonilization calculations 
the 3 and 4-bunches were confined to 12 sites on a ring, 
while for the 2-bunch we used 200 sites. This is sufficient 
to get accurate answers because the bunches are tightly 
bound with a mean particle spacing between 1 and 1.5. 
We have also computed the ground state energies us- 
ing the Green's Function Monte Carlo (GFMC) method 
, which allows exact calculations for many-boson sys- 
tems, for bunches confined on lattices with up to 100 
sites. Our GFMC calculations, which we describe below, 
verified that the size effect in the exact diagonalization 
was negligible || . 

We now consider the effects of including bunch-bunch 
interactions on the phase diagrams obtained in the EDL. 
Note that the contribution to the free energy per step due 
to bunch interactions is smaller than the leading contri- 
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bution /(rib, T) by a factor s 2 <C 1 at all temperatures. 
However, in the region around the phase boundaries in 
Fig. [j], called the "critical region", the difference in the 
step free energy in the rib and rib — 1 phases is small. 
When the bunch-bunch interaction energy becomes com- 
parable to this difference, we have to explicitly include 
it in the analysis in a non-perturbative way. Far from 
this region the bunch interactions can be treated per- 
turbatively. When the bunch interactions are taken into 
account, we anticipate on general grounds that the de- 
generacy of the step free energies in the rif, and rib — 1 
phases at the transition point is lifted, causing the bunch- 
ing transitions to become continuous. The width of the 
critical region, denoted by AT cr n oc s p ™'> , vanishes in 
the limit of vanishing densities (the EDL) as the bunch 
interactions approach zero in this limit. The exponent 
p nb , along with the shift in bunching transition temper- 
ature from the EDL due to bunch interactions, will be 
determined below. 

The bunch-bunch interaction energy in regions away 
from the critical region can be computed perturbatively 
by noting two key points: (l)The average spacing in a 
bunch (for 2-, 3- and 4-bunch phases in Fig.Jj], for u/g < 3, 
remains well under two lattice spacings M, so that the 
steps in a bunch remain tightly bound as a single entity. 
Since the bunches are very strongly bound, the stiffness of 
a bunch can then be approximated as = rif,E, where S 
is the stiffness of a single step. This expression for bunch 
stiffness is consistent with the experimental observations 
of the stiffnesses of the double, triple and quadruple steps 
made by Sudoh et al 0]. (2) The bunches interact with 
each other with a renormalized inverse square potential of 
strength gb = nig. This approximation is justified when 
the spacing between the bunches is very large. We note 
that the contributions to the energy of the system aris- 
ing from the short range bunch-bunch interactions are 
smaller than the inverse square contribution by a fac- 
tor proportional to the density of steps s <C 1. This 
can be shown very easily using a Hartree-Fock estima- 
tion of the short-ranged contribution M. The energy of 
the system can then be computed using the Calagero- 
Suthcrland model of hard-core bosons interacting via an 
inverse square law, for which an exact solution was pro- 
vided by Sutherland ||XC(] . Using his solution, we can write 
the energy per site in a system with bunches of size rib 
as 



E [s,T,n b ] _ f(nb,T)a x 
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| 1 + y/l + 2n|g . Using this energy, one 
can compute such quantities as surface stiffness and crys- 
tal shapes. These computations will appear elsewhere. |J 
In the critical region, since the energies of both the rib 



and rib — 1 bunches are nearly equal, we retain the terms 
of order s 3 arising from bunch-bunch interactions. In or- 
der to estimate the width of the critical region, we use 
the criterion Eo(s,T,rib — 1) ~ Eo(s,T,rib), to see that 
s 2 ~ \ f( n b, T)—f(rib—l, T)\ holds in the critical region. If 
|/K T) - f(n b - 1, T)| ~ |T- T c (n b -> n b - 1; 0)\ a ^ , we 
obtain the result p Ub = 2/a rib , for p nb introduced earlier. 
For bunching transitions first order in the EDL, we have 
a nb = 1 and therefore p nb = 2. As pointed out earlier, 
since all the bunching transitions shown in Fig. |l| were 
found to be first order in the EDL, the width of the crit- 
ical region is AT cni cx s 2 . Extending the phase diagram 
to larger values of u/g, we find that for u/g > 2.4 M, the 

1- bunch to 3-bunch transition preempts the 1-bunch to 

2- bunch transition and that for large u/g all transitions 
asymptotically approach the line l/g — 2(u/g). Further, 
all the bunching transitions that the system undergoes 
for u/g up to 8 were found to be first order Note 
next that the simplest conjecture for the shift of the tran- 
sition temperature due to bunch interactions is that it is 
of the same order of magnitude as the width of the crit- 
ical region. Explicitly, for bunching transitions from the 
rib-bunch to the n b -bunch phase that are first order in 
the EDL, we have 
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We have also performed GFMC simulations to study 
the physics of bunching transitions in the "thermody- 
namic limit" where the bunch interactions are treated 
non-perturbatively and to verify the results on bunching 
from exact diagonalization. To explicitly illustrate the 
bunching transitions, we compute the mixed estimator 
H of the pair correlation function, which for N bosons 
on L sites (L even and I < x < L/2 — 1) defined by 



(5) 



where x.i and the coordinates of the bosons i and 

i' and the brackets indicate the quantum average in the 
ground state. In Fig. ||we show plots of g[x) for u/g = 2. 
The figure vividly demonstrates the transitions from the 
1-bunch to 2-bunch, the 2-bunch to 3-bunch and the 3- 
bunch to 4-bunch. For comparison we have marked in the 
phase diagram in Fig. [j]the bunch sizes for each value of g 
obtained by counting the number of peaks in g(x). The 
phases predicted for 12 bosons confined on 120 sites is 
in excellent agreement with the EDL phase diagram ob- 
tained using exact diagonalization. Note that for points 
very close to the phase boundaries {g = 4 near the 3-4 
boundary and g = 2.96 near the 2-3 boundary) we see 
a continuous evolution from the (rif, — l)-bunch phase to 
the rib-bunch phase, consistent with lifting of level de- 
generacy as a result of the bunch-bunch interactions. We 
have studied the pair correlation functions for a number 
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of other values of u/g to verify the accuracy of the phase 
diagram shown in Fig. | [§. Our GFMC energies agree 
with Eq. (§). 
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FIG. 2. Boson pair correlation function g{x) (see Eq. (|5|)) 
for 12 bosons on a ring with 120 sites. We have chosen 
u/g — 2.0 and plotted g(x) for values of g indicated in the fig- 
ure. The points are also marked in Fig. |l| for comparison. By 
counting the number of peaks in g(x) and making use of the 
periodic boundary condition the number of bosons in a bunch 
can be established. We see that there are 3 4-bunches for 
g = 7.0, 4 3-bunches for g = 3.5 and 6 2-bunches for g = 2.85. 
For g = 2.25 we have a completely unbound 1-bunch phase, 
while the points g = 4 and g — 2.96 are in critical regions. 
The statistical error bars have only been shown for the case 
g = 7; they have similar magnitudes in all the other cases. 

Recently the effect of competing short and long range 
interactions on the finite temperature phase transitions 
on vicinal surfaces was considered by Lassig || and Bhat- 
tacharjee who used a continuum version of the model 
employed here. They analyzed the phase transition using 
renormalization group (RG) techniques, finding that, for 
a given g, when the short range attraction is sufficiently 
attractive, the steps tend to phase separate. This picture 
is fundamentally different from ours. In our picture, the 
steps do not phase separate. Instead, they form sepa- 
rated bunches of a finite size. The bunches then interact 
with each other with renormalized interaction parame- 
ters as discussed above. With increasing temperature 
the steps undergo a series of peeling transitions where the 
bunch size progressively decreases. The size of a bunch 
and peeling transitions are crucial effects of a physical 
region where the continuum model is not valid. 

In the experiments of Song et al M, above a tem- 
perature T c (s), the surface is uniformly stepped, while 
below T c (s) it consists of facets separated by bunched 
steps. The number of steps in a bunch was seen to be 
about 22. From the measurements it was inferred that 
s ~ \T c (0)-T c (s)\P, with/3 = 0.42 ±.1. Our picture sug- 
gests that we identify the experimentally observed tran- 
sition as a 1-bunch to 22-bunch transition. With our 
simple choice of potential, this particular transition does 



not appear. However, we do observe, for u/g > 2.4 [Bj, 
that the 1-3 transition preempts the 1-2 transition. For 
a potential with next-nearest neighbor attractive inter- 
actions of sufficient strength, even at T = 0, a 1-bunch 
directly goes to a 3-bunch, bypassing the 2-bunch. For 
potentials with more complicated features, it can be ex- 
pected that a 1-22 transition will occur. For quantitative 
progress here, a more detailed understanding of the step- 
step interaction on Si(113) surfaces is required. 

From Eq. (Q), if a 1-22 transition is first order in the 
EDL, i.e., the free energy curves /(22,T) and /(1,T) 
intersect with a finite slope, the shift in transition tem- 
perature scales like s 2 . It then follows that (3 = 0.5, in 
agreement with the experimentally observed exponent. 
It is important to note that in our picture describes 
a curve of continuous 1-22 bunch transitions rather than 
a curve of first order transitions associated with a tri- 
critical point. Within the recent RG calculations, (3 is 
non-universal, requiring a particular choice of u/g to pro- 
duce agreement with experiment. In contrast, our result 
[3 = 0.5 is robust, applying to any transition that is first- 
order in the EDL. All transitions in our model are indeed 
first-order in the EDL. It is precisely in this situation that 
the RG calculations fail. 

It remains difficult to reconcile some observed features 
of the bunch period and misorientation angle with our 
picture. In particular, the minimum observed in the 
bunch period as a function of the misorientation angle is 
readily understood jy within a picture, due to Marchenko 
p2| , where coexistence of two macroscopic facets is un- 
stable and replaced by a periodic groove structure which 
is in fact a bunch phase. The difference between the 
Marchenko picture and ours is that the bunch size in the 
former is determined by a balance between edge ener- 
gies and strain energies, whereas in the latter the bunch 
size is determined by the competition between short- and 
long-range step interactions. Further research, extend- 
ing our work to the case of large bunches not necessarily 
widely separated, is required to elucidate the connection 
between these two pictures. 

Our interpretation of the results of Song et al gains 
considerable support by the direct observation of single, 
double, triple, and quadruple steps by Sudoh et al. The 
bunch sizes in the two experiments are expected to dif- 
fer because the vicinal surfaces involved have different 
orientations with respect to the [113] direction, leading 
to different short and long range step-step interactions. 
Further exploration of the vicinal surfaces of Si near the 
[113] direction in an effort to fully elucidate the connec- 
tions between the two sets of experiments will be crucial 
to the further development of this area. 
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